UTS MPDW

Angga Fathan Rofiqy

06 October, 2023

Kode di Hide dalam default, untuk menampilkan kode, klik Code .

#                      -=( Install & Load Package Function )=-
install_load <- function (package1, ...)  {   

   # convert arguments to vector
   packages <- c(package1, ...)

   # start loop to determine if each package is installed
   for(package in packages){

       # if package is installed locally, load
       if(package %in% rownames(installed.packages()))
          do.call('library', list(package))

       # if package is not installed locally, download, then load
       else {
          install.packages(package)
          do.call("library", list(package))
       }
   } 
}

path <- function(){
  gsub  ( "\\\\",  "/",  readClipboard ()  )
}
#Copy path, Panggil function di console
#Copy r path, paste ke var yang diinginkan

0.1 Data Preparation

0.1.1 Import Data

install_load('rio')
raw.data <- import("https://raw.githubusercontent.com/Zen-Rofiqy/STA1341-MPDW/main/UTS/Data%20UTS%20Fix.csv")
data <- raw.data

0.2 Packages

install_load('ggplot2','tsibble','tsibble','MASS','forecast','TSA','TTR','aTSA','graphics')

## Warning: package 'ggplot2' was built under R version 4.2.3

## Warning: package 'tsibble' was built under R version 4.2.3

## 
## Attaching package: 'tsibble'

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, union

## Warning: package 'forecast' was built under R version 4.2.3

## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

## Warning: package 'TSA' was built under R version 4.2.3

## Registered S3 methods overwritten by 'TSA':
##   method       from    
##   fitted.Arima forecast
##   plot.Arima   forecast

## 
## Attaching package: 'TSA'

## The following objects are masked from 'package:stats':
## 
##     acf, arima

## The following object is masked from 'package:utils':
## 
##     tar

## 
## Attaching package: 'aTSA'

## The following object is masked from 'package:forecast':
## 
##     forecast

## The following object is masked from 'package:graphics':
## 
##     identify

0.3 Data

set.seed(575)
ma2 <- arima.sim(list(order = c(0,0,2), ma = c(0.55,0.65)), n = 175)

Data kemudian dibagi menjadi data latih dan data uji. Pembagian kali ini dilakukan dengan proporsi / perbandingan, yaitu 80:20.

data.train <- data[1: round(nrow(data) *80/100),]
data.test <- data[round(nrow(data) *80/100): nrow(data),]
data.train.ts <- ts(data.train[,1])
data.test.ts <- ts(data.test[,1])

0.3.1 Eksplorasi Data

Sebelum masuk dalam tahap pemodelan, dilakukan eksplorasi data dengan plot deret waktu untuk melihat pola data.

data.ts <- raw.data[,2]
ts.plot(data.ts, xlab="Time Period", ylab="Data", 
        main = "Time Series Data", col='orange', lwd=2)
points(data.ts, col='orange', lwd=1.5)

ctsa.train_test <-
ggplot() + 
  geom_line(data = data.train, 
            aes(x = Periode, y = Data, col = "Data Latih")) +
  geom_line(data = data.test, 
            aes(x = Periode, y = Data, col = "Data Uji")) +
  labs(x = "Periode Waktu", y = "Data", color = "Legend") +
  scale_colour_manual(name="Keterangan:", 
                      breaks = c("Data Latih", "Data Uji"),
                      values = c("orange", "cyan4")) + 
  theme_bw() + theme(legend.position = "bottom",
                     plot.caption = element_text(hjust=0.5, size=12))
ctsa.train_test

## Data Bangkitan

0.3.2 Pembangkitan Data

Data yang akan dibangkitkan adalah data dengan model MA(2) sebagai berikut.

set.seed(575)
ma2 <- arima.sim(list(order = c(0,0,2), ma = c(0.55,0.65)), n = 175)

Data kemudian dibagi menjadi data latih dan data uji. Pembagian kali ini dilakukan dengan proporsi / perbandingan, yaitu 80:20.

ma2 <- ma2[-c(1:25)]
ma2.train <- ma2[1:120]
ma2.test <- ma2[121:150]

0.3.3 Eksplorasi Data

Sebelum masuk dalam tahap pemodelan, dilakukan eksplorasi data dengan plot deret waktu untuk melihat pola data.

#--PLOT TIME SERIES--#
plot(ma2.train,
     col = "navyblue",
     lwd = 1,
     type = "o",
     xlab = "Time",
     ylab = "Data")

Berdasarkan plot data deret waktu di atas, terlihat data cenderung stasioner dalam rataan dan ragam. Data stasioner dalam rataan karena menyebar/bergerak di sekitar nilai tengahnya (0) dan dikatakan stasioner dalam ragam karena memiliki lebar pita yang cenderung sama. Selain dengan plot data deret waktu, akan dilakukan pengecekan stasioneritas data dengan plot ACF dan uji ADF.

#--CEK KESTASIONERAN---#
acf(ma2.train, main="ACF", lag.max=20)

Berdasarkan plot ACF di atas, dapat dilihat bahwa plot cuts off pada lag ke-2. Hal ini sesuai dengan proses pembangkitan model MA(2).

adf.test(ma2.train)

## Augmented Dickey-Fuller Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##      lag   ADF p.value
## [1,]   0 -7.82    0.01
## [2,]   1 -5.12    0.01
## [3,]   2 -6.26    0.01
## [4,]   3 -4.62    0.01
## [5,]   4 -3.46    0.01
## Type 2: with drift no trend 
##      lag   ADF p.value
## [1,]   0 -7.82  0.0100
## [2,]   1 -5.11  0.0100
## [3,]   2 -6.26  0.0100
## [4,]   3 -4.63  0.0100
## [5,]   4 -3.46  0.0116
## Type 3: with drift and trend 
##      lag   ADF p.value
## [1,]   0 -8.09  0.0100
## [2,]   1 -5.50  0.0100
## [3,]   2 -6.77  0.0100
## [4,]   3 -5.12  0.0100
## [5,]   4 -4.01  0.0115
## ---- 
## Note: in fact, p.value = 0.01 means p.value <= 0.01

#stasioner

\(H_0\) : Data tidak stasioner dalam rataan

\(H_1\) : Data stasioner dalam rataan

Berdasarkan uji ADF tersebut, didapat p-value sebesar 0.01358 yang lebih kecil dari taraf nyata 5% sehingga tolak \(H_0\) dan menandakan bahwa data stasioner dalam rataan. Hal ini sesuai dengan hasil eksplorasi menggunakan plot time series dan plot ACF.

0.3.4 Spesifikasi Model

#---SPESIFIKASI MODEL---#
par(mfrow = c(1,2))
acf(ma2.train, main="ACF", lag.max=20) #ARIMA(0,0,2)
pacf(ma2.train, main="PACF", lag.max=20) #ARIMA(1,0,0)

par(mfrow = c(1,1))

Berdasarkan Plot ACF, terlihat cuts off pada lag ke-2 sehingga dapat kita asumsikan model yang terbentuk adalah ARIMA(0,0,2). Selanjutnya, berdasarkan plot PACF, terlihat cuts off pada lag pertama sehingga model yang terbentuk adalah ARIMA(1,0,0). Selain dengan plot ACF dan PACF, penentuan spesifikasi model dilakukan dengan extended ACF (EACF) berikut ini.

eacf(ma2.train)

## AR/MA
##   0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 x x o o o o o o o o o  o  o  o 
## 1 x x o o o o o o o o o  o  o  o 
## 2 x x o o o o o o o o o  o  o  o 
## 3 x x x x o o o o o o o  o  o  o 
## 4 x x o o o o o o o o o  o  o  o 
## 5 x x o o o o o o o o o  o  o  o 
## 6 x x x o o o o o o o o  o  o  o 
## 7 x x o o o o o o o o o  o  o  o

#ARIMA(0,0,2) #ARIMA(1,0,3) #ARIMA(2,0,3) #ARIMA(3,0,3)
#Terdapat 5 model tentatif

Menggunakan plot EACF, dapat diambil beberapa model dengan melihat ujung segitiga yang terbentuk, antara lain ARIMA(0,0,2), ARIMA(1,0,3), ARIMA(2,0,3), dan ARIMA(3,0,3).

0.3.5 Pendugaan Parameter

Selanjutnya akan dilakukan pendugaan parameter kelima model ARIMA yang terbentuk sebelumnya. Pendugaan dilakukan dengan fungsi Arima() yang dilanjutkan dengan melihat nilai AIC pada ringkasan data dan melihat signifikansi parameter.

#---PENDUGAAN PARAMETER MODEL---#
model1.ma2=Arima(ma2.train, order=c(0,0,2),method="ML")
summary(model1.ma2) #AIC=326.87

## Series: ma2.train 
## ARIMA(0,0,2) with non-zero mean 
## 
## Coefficients:
##          ma1     ma2    mean
##       0.5042  0.6656  0.1198
## s.e.  0.0838  0.0802  0.2105
## 
## sigma^2 = 1.175:  log likelihood = -179.04
## AIC=366.08   AICc=366.43   BIC=377.23
## 
## Training set error measures:
##                        ME     RMSE       MAE      MPE     MAPE      MASE
## Training set -0.004465907 1.070145 0.8639788 154.2116 248.8555 0.7108001
##                    ACF1
## Training set -0.1185747

lmtest::coeftest(model1.ma2) #seluruh parameter signifikan

## 
## z test of coefficients:
## 
##           Estimate Std. Error z value  Pr(>|z|)    
## ma1       0.504177   0.083802  6.0163 1.785e-09 ***
## ma2       0.665552   0.080152  8.3036 < 2.2e-16 ***
## intercept 0.119801   0.210531  0.5690    0.5693    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

model2.ma2=Arima(ma2.train, order=c(1,0,0),method="ML")
summary(model2.ma2) #AIC=340.47

## Series: ma2.train 
## ARIMA(1,0,0) with non-zero mean 
## 
## Coefficients:
##          ar1    mean
##       0.3179  0.1224
## s.e.  0.0864  0.1633
## 
## sigma^2 = 1.526:  log likelihood = -194.67
## AIC=395.33   AICc=395.54   BIC=403.7
## 
## Training set error measures:
##                        ME     RMSE       MAE      MPE     MAPE      MASE
## Training set -0.004054871 1.224874 0.9882004 101.8339 166.6078 0.8129979
##                     ACF1
## Training set -0.07888664

lmtest::coeftest(model2.ma2) #seluruh parameter signifikan

## 
## z test of coefficients:
## 
##           Estimate Std. Error z value  Pr(>|z|)    
## ar1        0.31794    0.08644  3.6782 0.0002349 ***
## intercept  0.12241    0.16331  0.7496 0.4535097    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

model3.ma2=Arima(ma2.train, order=c(1,0,3),method="ML")
summary(model3.ma2) #AIC=329.22

## Series: ma2.train 
## ARIMA(1,0,3) with non-zero mean 
## 
## Coefficients:
##          ar1     ma1     ma2      ma3    mean
##       0.0199  0.3647  0.5710  -0.2140  0.1180
## s.e.  0.4065  0.3948  0.2188   0.2809  0.1676
## 
## sigma^2 = 1.149:  log likelihood = -176.79
## AIC=365.58   AICc=366.33   BIC=382.31
## 
## Training set error measures:
##                        ME     RMSE       MAE      MPE     MAPE     MASE
## Training set -0.003482706 1.049274 0.8648842 173.1621 259.3462 0.711545
##                       ACF1
## Training set -0.0006138296

lmtest::coeftest(model3.ma2) #tidak ada yang signifikan

## 
## z test of coefficients:
## 
##            Estimate Std. Error z value Pr(>|z|)   
## ar1        0.019876   0.406479  0.0489 0.961001   
## ma1        0.364683   0.394772  0.9238 0.355601   
## ma2        0.570966   0.218840  2.6091 0.009079 **
## ma3       -0.214042   0.280904 -0.7620 0.446074   
## intercept  0.117990   0.167626  0.7039 0.481501   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

model4.ma2=Arima(ma2.train, order=c(2,0,3),method="ML")
summary(model4.ma2) #AIC=330.6

## Series: ma2.train 
## ARIMA(2,0,3) with non-zero mean 
## 
## Coefficients:
##           ar1      ar2     ma1     ma2      ma3    mean
##       -0.0752  -0.0332  0.4597  0.6374  -0.1450  0.1183
## s.e.   1.0018   0.3061  0.9988  0.6720   0.7244  0.1680
## 
## sigma^2 = 1.159:  log likelihood = -176.79
## AIC=367.57   AICc=368.57   BIC=387.08
## 
## Training set error measures:
##                        ME     RMSE       MAE      MPE     MAPE     MASE
## Training set -0.003608762 1.049191 0.8643288 173.1499 257.9485 0.711088
##                       ACF1
## Training set -0.0008176693

lmtest::coeftest(model4.ma2) #hanya ma2 yang signifikan

## 
## z test of coefficients:
## 
##           Estimate Std. Error z value Pr(>|z|)
## ar1       -0.07515    1.00181 -0.0750   0.9402
## ar2       -0.03320    0.30611 -0.1085   0.9136
## ma1        0.45973    0.99883  0.4603   0.6453
## ma2        0.63743    0.67200  0.9486   0.3428
## ma3       -0.14504    0.72436 -0.2002   0.8413
## intercept  0.11826    0.16804  0.7037   0.4816

model5.ma2=Arima(ma2.train, order=c(3,0,3),method="ML")
summary(model5.ma2) #AIC=329.87

## Series: ma2.train 
## ARIMA(3,0,3) with non-zero mean 
## 
## Coefficients:
##           ar1      ar2      ar3     ma1     ma2     ma3    mean
##       -1.1458  -0.3249  -0.1172  1.5357  1.3411  0.6916  0.1165
## s.e.   0.2499   0.2039   0.1281  0.2359  0.2118  0.1799  0.1680
## 
## sigma^2 = 1.164:  log likelihood = -176.52
## AIC=369.05   AICc=370.35   BIC=391.35
## 
## Training set error measures:
##                        ME     RMSE       MAE      MPE    MAPE     MASE
## Training set -0.002222695 1.046847 0.8608718 163.7022 258.074 0.708244
##                      ACF1
## Training set -0.002236272

lmtest::coeftest(model5.ma2) #hanya ma1 dan ma2 yang signifikan

## 
## z test of coefficients:
## 
##           Estimate Std. Error z value  Pr(>|z|)    
## ar1       -1.14580    0.24987 -4.5856 4.526e-06 ***
## ar2       -0.32493    0.20394 -1.5933 0.1110965    
## ar3       -0.11718    0.12809 -0.9148 0.3602828    
## ma1        1.53566    0.23593  6.5090 7.566e-11 ***
## ma2        1.34111    0.21183  6.3309 2.437e-10 ***
## ma3        0.69157    0.17987  3.8449 0.0001206 ***
## intercept  0.11646    0.16797  0.6933 0.4881021    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#model yang dipilih adalah model 1, yaitu ARIMA(0,0,2)

Berdasarkan pendugaan parameter di atas, nilai AIC terkecil dimiliki oleh model ARIMA(0,0,2) dan parameter model ARIMA(0,0,2) juga seluruhnya signifikan sehingga model yang dipilih adalah model ARIMA(0,0,2).

0.3.6 Analisis Sisaan

Model terbaik hasil identifikasi kemudian dicek asumsi sisaannya. Sisaan model ARIMA harus memenuhi asumsi normalitas, kebebasan, dan kehomogenan ragam. Diagnostik model dilakukan secara eksplorasi dan uji formal.

0.3.6.1 Eksplorasi Sisaan

#Eksplorasi
sisaan.ma2 <- model1.ma2$residuals
par(mfrow=c(2,2))
qqnorm(sisaan.ma2)
qqline(sisaan.ma2, col = "blue", lwd = 2)
plot(c(1:length(sisaan.ma2)),sisaan.ma2)
acf(sisaan.ma2)
pacf(sisaan.ma2)

par(mfrow = c(1,1))

Berdasarkan plot kuantil-kuantil normal, secara eksplorasi ditunjukkan sisaan menyebar normal mengikuti garis \(45^{\circ}\). Kemudian dapat dilihat juga lebar pita sisaan yang cenderung sama menandakan bahwa sisaan memiliki ragam yang homogen. Akan tetapi, plot ACF dan PACF sisaan ARIMA(0,0,2) signifikan pada lag ke-6 sehingga sisaan tidak saling bebas. Kondisi ini akan diuji lebih lanjut dengan uji formal.

0.3.6.2 Uji Formal

#1) Sisaan Menyebar Normal
ks.test(sisaan.ma2,"pnorm")

## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  sisaan.ma2
## D = 0.067078, p-value = 0.6528
## alternative hypothesis: two-sided

#tak tolak H0 > sisaan menyebar normal

Selain dengan eksplorasi, asumsi tersebut dapat diuji menggunakan uji formal. Pada tahapan ini uji formal yang digunakan untuk normalitas adalah uji Kolmogorov-Smirnov (KS). Hipotesis pada uji KS adalah sebagai berikut.

\(H_0\) : Sisaan menyebar normal

\(H_1\) : Sisaan tidak menyebar normal

Berdasarkan uji KS tersebut, didapat p-value sebesar 0.9788 yang lebih besar dari taraf nyata 5% sehingga tak tolak \(H_0\) dan menandakan bahwa sisaan menyebar normal. Hal ini sesuai dengan hasil eksplorasi menggunakan plot kuantil-kuantil normal.

#2) Sisaan saling bebas/tidak ada autokorelasi
Box.test(sisaan.ma2, type = "Ljung")

## 
##  Box-Ljung test
## 
## data:  sisaan.ma2
## X-squared = 1.7297, df = 1, p-value = 0.1884

#tak tolak H0 > sisaan saling bebas

Selanjutnya akan dilakukan uji formal untuk kebebasan sisaan menggunakan uji Ljung-Box. Hipotesis yang digunakan adalah sebagai berikut.

\(H_0\) : Sisaan saling bebas

\(H_1\) : Sisaan tidak tidak saling bebas

Berdasarkan uji Ljung-Box tersebut, didapat p-value sebesar 0.5082 yang lebih besar dari taraf nyata 5% sehingga tak tolak \(H_0\) dan menandakan bahwa sisaan saling bebas. Hal ini berbeda dengan eksplorasi.

#3) Sisaan homogen
Box.test((sisaan.ma2)^2, type = "Ljung")

## 
##  Box-Ljung test
## 
## data:  (sisaan.ma2)^2
## X-squared = 0.18307, df = 1, p-value = 0.6687

#tak tolak H0 > sisaan homogen

Hipotesis yang digunakan untuk uji kehomogenan ragam adalah sebagai berikut.

\(H_0\) : Ragam sisaan homogen

\(H_1\) : Ragam sisaan tidak homogen

Berdasarkan uji Ljung-Box terhadap sisaan kuadrat tersebut, didapat p-value sebesar 0.116 yang lebih besar dari taraf nyata 5% sehingga tak tolak \(H_0\) dan menandakan bahwa ragam sisaan homogen.

#4) Nilai tengah sisaan sama dengan nol
t.test(sisaan.ma2, mu = 0, conf.level = 0.95)

## 
##  One Sample t-test
## 
## data:  sisaan.ma2
## t = -0.045524, df = 119, p-value = 0.9638
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -0.1987121  0.1897803
## sample estimates:
##    mean of x 
## -0.004465907

#tak tolak h0 > nilai tengah sisaan sama dengan 0

Terakhir, dengan uji-t, akan dicek apakah nilai tengah sisaan sama dengan nol. Hipotesis yang diujikan sebagai berikut.

\(H_0\) : nilai tengah sisaan sama dengan 0

\(H_1\) : nilai tengah sisaan tidak sama dengan 0

Berdasarkan uji-ttersebut, didapat p-value sebesar 0.9594 yang lebih besar dari taraf nyata 5% sehingga tak tolak \(H_0\) dan menandakan bahwa nilai tengah sisaan sama dengan nol. Hal ini berbeda dengan eksplorasi.

0.3.7 Overfitting

Tahapan selanjutnya adalah overfitting dilakukan dengan menaikkan orde AR(p) dan MA(q) dari model ARIMA(0,0,2) untuk melihat apakah terdapat model lain yang lebih baik dari model saat ini. Kandidat model overfitting adalah ARIMA(1,0,2) dan ARIMA(0,0,3).

#---OVERFITTING---#
model1a.ma2=Arima(ma2.train, order=c(1,0,2),method="ML")
summary(model1a.ma2) #327.31

## Series: ma2.train 
## ARIMA(1,0,2) with non-zero mean 
## 
## Coefficients:
##           ar1     ma1     ma2    mean
##       -0.2545  0.6491  0.7239  0.1207
## s.e.   0.1179  0.0809  0.0753  0.1805
## 
## sigma^2 = 1.143:  log likelihood = -177
## AIC=364   AICc=364.53   BIC=377.94
## 
## Training set error measures:
##                        ME     RMSE       MAE      MPE     MAPE     MASE
## Training set -0.004810005 1.051299 0.8679072 168.4309 258.7567 0.714032
##                     ACF1
## Training set -0.01048925

lmtest::coeftest(model1a.ma2) #ar1 tidak signifikan

## 
## z test of coefficients:
## 
##            Estimate Std. Error z value  Pr(>|z|)    
## ar1       -0.254474   0.117933 -2.1578   0.03095 *  
## ma1        0.649118   0.080942  8.0195 1.061e-15 ***
## ma2        0.723925   0.075342  9.6085 < 2.2e-16 ***
## intercept  0.120712   0.180540  0.6686   0.50374    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

model1b.ma2=Arima(ma2.train, order=c(0,0,3),method="ML")
summary(model1b.ma2) #327.24

## Series: ma2.train 
## ARIMA(0,0,3) with non-zero mean 
## 
## Coefficients:
##          ma1     ma2      ma3    mean
##       0.3836  0.5810  -0.2010  0.1182
## s.e.  0.0888  0.0807   0.0924  0.1682
## 
## sigma^2 = 1.139:  log likelihood = -176.79
## AIC=363.59   AICc=364.11   BIC=377.52
## 
## Training set error measures:
##                        ME     RMSE       MAE      MPE     MAPE      MASE
## Training set -0.003598688 1.049285 0.8651174 172.9022 259.0877 0.7117368
##                     ACF1
## Training set 0.000279606

lmtest::coeftest(model1b.ma2) #ma3 tidak signifikan

## 
## z test of coefficients:
## 
##            Estimate Std. Error z value  Pr(>|z|)    
## ma1        0.383565   0.088791  4.3199 1.561e-05 ***
## ma2        0.580966   0.080676  7.2012 5.966e-13 ***
## ma3       -0.200952   0.092431 -2.1741    0.0297 *  
## intercept  0.118227   0.168194  0.7029    0.4821    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#model yang dipilih adalah model awal, yaitu ARIMA(0,0,2)

Berdasarkan kedua model hasil overfitting di atas, model ARIMA(1,0,2) dan ARIMA(0,0,3) memiliki AIC yang lebih besar dibandingkan dengan model ARIMA(0,0,2) dan parameter kedua model ARIMA(1,0,2) dan ARIMA(0,0,3) tidak seluruhnya signifikan. Oleh karena itu, model ARIMA(0,0,2) akan tetap digunakan untuk melakukan peramalan.

0.3.8 Peramalan

Peramalan dilakukan menggunakan fungsi forecast() . Contoh peramalan berikut ini dilakukan untuk 30 hari ke depan.

#---FORECAST---#
ramalan <- forecast::forecast(model1.ma2, h = 30)
ramalan

##     Point Forecast     Lo 80    Hi 80     Lo 95    Hi 95
## 121    -0.01625282 -1.405170 1.372665 -2.140419 2.107913
## 122    -0.49156069 -2.047021 1.063900 -2.870432 1.887310
## 123     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 124     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 125     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 126     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 127     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 128     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 129     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 130     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 131     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 132     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 133     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 134     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 135     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 136     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 137     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 138     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 139     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 140     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 141     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 142     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 143     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 144     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 145     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 146     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 147     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 148     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 149     0.11980094 -1.689610 1.929212 -2.647454 2.887056
## 150     0.11980094 -1.689610 1.929212 -2.647454 2.887056

data.ramalan <- ramalan$mean
plot(ramalan)

Berdasarkan hasil plot ramalan di atas, dapat dilihat bahwa ramalan ARIMA(0,0,2) cenderung meningkat di awal periode dan stabil hingga akhir periode. Selanjutnya, dapat dicari nilai akurasi antara hasil ramalan dengan data uji sebagai berikut.

perbandingan<-matrix(data=c(ma2.test, data.ramalan),
                     nrow = 30, ncol = 2)
colnames(perbandingan)<-c("Aktual","Hasil Forecast")
perbandingan

##            Aktual Hasil Forecast
##  [1,]  1.10823525    -0.01625282
##  [2,] -0.19890485    -0.49156069
##  [3,] -0.69806908     0.11980094
##  [4,] -1.86358397     0.11980094
##  [5,] -1.66113761     0.11980094
##  [6,] -1.28449426     0.11980094
##  [7,] -0.58315721     0.11980094
##  [8,]  0.31073371     0.11980094
##  [9,]  0.28744684     0.11980094
## [10,]  0.11804291     0.11980094
## [11,] -0.37690166     0.11980094
## [12,]  0.61037938     0.11980094
## [13,]  0.55668732     0.11980094
## [14,]  1.83166574     0.11980094
## [15,]  1.35978016     0.11980094
## [16,]  0.35996629     0.11980094
## [17,]  0.24773649     0.11980094
## [18,] -1.98310931     0.11980094
## [19,] -0.89190922     0.11980094
## [20,] -2.02462962     0.11980094
## [21,] -1.17509140     0.11980094
## [22,] -1.21066999     0.11980094
## [23,] -2.06176773     0.11980094
## [24,] -3.82058583     0.11980094
## [25,] -2.22554173     0.11980094
## [26,] -1.68222822     0.11980094
## [27,]  0.07758707     0.11980094
## [28,] -1.02825180     0.11980094
## [29,] -0.15833731     0.11980094
## [30,] -2.89955052     0.11980094

accuracy(data.ramalan, ma2.test)

##                  ME     RMSE      MAE      MPE     MAPE
## Test set -0.7935424 1.526591 1.195085 84.36333 97.89875

0.4 Data Asli

Digunakan data kurs yang dalam hal ini hanya digunakan data 500 periode awal

datakurs<-raw.data[,2]
datakurs.ts<-ts(datakurs)

0.4.1 Eksplorasi Data

0.4.1.1 Plot Data Penuh

plot.ts(datakurs.ts, lty=1, xlab="waktu", ylab="Kurs", main="Plot Data Kurs")

Berdasarkan plot data deret waktu, terlihat bahwa data cenderung memiliki trend yang naik. Berdasarkan pola data, pembagian data latih dan data uji ditetapkan dengan proporsi 86%:14%.

0.4.1.2 Plot Data Latih

kurstrain<-datakurs[1:120]
train.ts<-ts(kurstrain)
plot.ts(train.ts, lty=1, xlab="waktu", ylab="Kurs", main="Plot Kurs Train")

Berdasarkan plot data deret waktu pada data latih, terlihat bahwa data cenderung memiliki trend yang naik dan cenderung tidak bergerak pada nilai tengah tertentu. Hal ini mengindikasikan bahwa data tidak stasioner dalam rataan.

0.4.1.3 Plot Data Uji

kurstest<-datakurs[120:151]
test.ts<-ts(kurstest)
plot.ts(test.ts, lty=1, xlab="waktu", ylab="Kurs", main="Plot Kurs Test")

0.4.2 Uji Stasioneritas Data

0.4.2.1 Plot ACF

acf(train.ts)

Berdasarkan plot ACF, terlihat bahwa plot ACF data menurun secara perlahan (tails of slowly). Hal ini juga menjadi indikasi bahwa data tidak stasioner dalam rataan

0.4.2.2 Uji ADF

tseries::adf.test(train.ts)

## 
##  Augmented Dickey-Fuller Test
## 
## data:  train.ts
## Dickey-Fuller = -2.5599, Lag order = 4, p-value = 0.344
## alternative hypothesis: stationary

\(H_0\) : Data tidak stasioner dalam rataan

\(H_1\) : Data stasioner dalam rataan

Berdasarkan uji ADF tersebut, didapat p-value sebesar 0.5553 yang lebih besar dari taraf nyata 5% sehingga tak tolak \(H_0\) dan menandakan bahwa data tidak stasioner dalam rataan. Hal ini sesuai dengan hasil eksplorasi menggunakan plot time series dan plot ACF, sehingga ketidakstasioneran model kedepannya harus ditangani

0.4.2.3 Plot Box-Cox

index <- seq(1:120)
bc = boxcox(train.ts~index, lambda = seq(5,10,by=1))

#Nilai Rounded Lambda
lambda <- bc$x[which.max(bc$y)]
lambda

## [1] 5

#SK
bc$x[bc$y > max(bc$y) - 1/2 * qchisq(.95,1)]

## [1] 5.000000 5.050505 5.101010

Plot Boxcox menunjukkan nilai rounded value (\(\lambda\)) optimum sebesar 6,64 dan pada selang kepercayaan 95% nilai memiliki batas bawah 0,48 dan batas atas 5,27. Selang tersebut memuat nilai satu sehingga dapat dikatakan bahwa data bangkitan stasioner dalam ragam.

0.4.3 Penanganan Ketidakstasioneran Data

train.diff<-diff(train.ts,differences = 1)
plot.ts(train.diff, lty=1, xlab="waktu", ylab="Data Difference 1 Kurs", main="Plot Difference Kurs")

Berdasarkan plot data deret waktu, terlihat bahwa data sudah stasioner dalam rataan ditandai dengan data bergerak pada nilai tengah tertentu (tidak terdapat trend ataupun musiman pada data)

0.4.3.1 Plot ACF

acf(train.diff)

Berdasarkan plot tersebut, terlihat bahwa plot ACF cuts off pada lag ke 1. Hal ini menandakan data sudah stasioner dalam rataan dan ketidakstasioneran data telah berhasil tertangani.

0.4.3.2 Uji ADF

tseries::adf.test(train.diff)

## Warning in tseries::adf.test(train.diff): p-value smaller than printed p-value

## 
##  Augmented Dickey-Fuller Test
## 
## data:  train.diff
## Dickey-Fuller = -4.5242, Lag order = 4, p-value = 0.01
## alternative hypothesis: stationary

\(H_0\) : Data tidak stasioner dalam rataan

\(H_1\) : Data stasioner dalam rataan

Berdasarkan uji ADF tersebut, didapat p-value sebesar 0.01 yang lebih kecil dari taraf nyata 5% sehingga tolak \(H_0\) atau data stasioner dalam rataan. Hal ini sesuai dengan hasil eksplorasi menggunakan plot time series dan plot ACF, sehingga dalam hal ini ketidakstasioneran data sudah berhasil ditangani dan dapat dilanjutkan ke pemodelan

0.4.4 Identifikasi Model

0.4.4.1 Plot ACF

acf(train.diff)

Berdasarkan plot tersebut, terlihat bahwa plot ACF cenderung cuts off pada lag ke 1, sehingga jika plot PACF dianggap tails of, maka model tentatifnya adalah ARIMA(0,1,1).

0.4.4.2 Plot PACF

pacf(train.diff)

Berdasarkan plot tersebut, terlihat bahwa plot PACF cenderung cuts off pada lag ke 1, sehingga jika plot ACF dianggap tails of, maka model tentatifnya adalah ARIMA(1,1,0).

Jika baik plot ACF maupun plot PACF keduanya dianggap tails of, maka model yang terbentuk adalah ARIMA(1,1,1)

0.4.4.3 Plot EACF

eacf(train.diff)

## AR/MA
##   0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 x o o o o x o o o o o  o  o  x 
## 1 x o o o o o o o o o o  o  o  o 
## 2 x o o o o o o x o o o  o  o  o 
## 3 x x o o o o o o o o o  o  o  o 
## 4 x o o o o o o o o o o  o  o  o 
## 5 x o x x x o o o o o o  o  o  o 
## 6 o x x x x o o o o o o  o  o  o 
## 7 x x x x o o x o o o o  o  o  o

Identifikasi model menggunakan plot EACF dilakukan dengan melihat ujung segitiga pada pola segitiga nol. Dalam hal ini model tentatif yang terbentuk adalah ARIMA(0,1,2), ARIMA(1,1,2), ARIMA(2,1,2), dan ARIMA(3,1,2).

0.4.5 Pendugaan Parameter Model Tentatif

0.4.5.1 ARIMA(0,1,1)

model1.da=Arima(train.diff, order=c(0,1,1),method="ML")
summary(model1.da) #AIC=4753.18

## Series: train.diff 
## ARIMA(0,1,1) 
## 
## Coefficients:
##           ma1
##       -1.0000
## s.e.   0.0229
## 
## sigma^2 = 9080:  log likelihood = -707.04
## AIC=1418.07   AICc=1418.18   BIC=1423.61
## 
## Training set error measures:
##                     ME     RMSE      MAE MPE MAPE      MASE       ACF1
## Training set 0.5237226 94.48326 72.30748 Inf  Inf 0.6249841 -0.1848378

lmtest::coeftest(model1.da) #seluruh parameter signifikan

## 
## z test of coefficients:
## 
##      Estimate Std. Error z value  Pr(>|z|)    
## ma1 -1.000000   0.022854 -43.756 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

0.4.5.2 ARIMA(1,1,0)

model2.da=Arima(train.diff, order=c(1,1,0),method="ML")
summary(model2.da) #AIC=4917.41

## Series: train.diff 
## ARIMA(1,1,0) 
## 
## Coefficients:
##           ar1
##       -0.5662
## s.e.   0.0766
## 
## sigma^2 = 14696:  log likelihood = -733.25
## AIC=1470.5   AICc=1470.61   BIC=1476.04
## 
## Training set error measures:
##                     ME     RMSE     MAE MPE MAPE     MASE      ACF1
## Training set -1.435216 120.2047 89.1368 Inf  Inf 0.770447 -0.249072

lmtest::coeftest(model2.da) #seluruh parameter signifikan

## 
## z test of coefficients:
## 
##      Estimate Std. Error z value Pr(>|z|)    
## ar1 -0.566243   0.076567 -7.3954 1.41e-13 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

0.4.5.3 ARIMA(1,1,1)

model3.da=Arima(train.diff, order=c(1,1,1),method="ML")
summary(model3.da) #AIC=4761.39

## Series: train.diff 
## ARIMA(1,1,1) 
## 
## Coefficients:
##           ar1     ma1
##       -0.1886  -1.000
## s.e.   0.0910   0.024
## 
## sigma^2 = 8809:  log likelihood = -704.93
## AIC=1415.86   AICc=1416.07   BIC=1424.17
## 
## Training set error measures:
##                   ME     RMSE      MAE MPE MAPE      MASE        ACF1
## Training set 1.53219 92.66424 71.10624 Inf  Inf 0.6146013 -0.01374354

lmtest::coeftest(model3.da) #seluruh parameter signifikan

## 
## z test of coefficients:
## 
##      Estimate Std. Error  z value Pr(>|z|)    
## ar1 -0.188595   0.090952  -2.0736  0.03812 *  
## ma1 -0.999994   0.024011 -41.6468  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

0.4.5.4 ARIMA(0,1,2)

model4.da=Arima(train.diff, order=c(0,1,2),method="ML")
summary(model4.da) #AIC=4748.3

## Series: train.diff 
## ARIMA(0,1,2) 
## 
## Coefficients:
##           ma1     ma2
##       -1.2001  0.2001
## s.e.   0.0932  0.0901
## 
## sigma^2 = 8777:  log likelihood = -704.77
## AIC=1415.54   AICc=1415.75   BIC=1423.85
## 
## Training set error measures:
##                    ME    RMSE      MAE MPE MAPE      MASE          ACF1
## Training set 1.609717 92.4947 71.04626 Inf  Inf 0.6140828 -0.0008224433

lmtest::coeftest(model4.da) #seluruh parameter signifikan

## 
## z test of coefficients:
## 
##      Estimate Std. Error  z value Pr(>|z|)    
## ma1 -1.200084   0.093230 -12.8723  < 2e-16 ***
## ma2  0.200085   0.090093   2.2209  0.02636 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

0.4.5.5 ARIMA(1,1,2)

model5.da=Arima(train.diff, order=c(1,1,2),method="ML")
summary(model5.da) #AIC=4749.85

## Series: train.diff 
## ARIMA(1,1,2) 
## 
## Coefficients:
##          ar1      ma1     ma2
##       0.0077  -1.2076  0.2076
## s.e.  0.3606   0.3494  0.3485
## 
## sigma^2 = 8853:  log likelihood = -704.77
## AIC=1417.54   AICc=1417.89   BIC=1428.62
## 
## Training set error measures:
##                    ME     RMSE     MAE MPE MAPE      MASE          ACF1
## Training set 1.609121 92.49418 71.0488 Inf  Inf 0.6141047 -0.0009336968

lmtest::coeftest(model5.da) #terdapat parameter tidak signifikan

## 
## z test of coefficients:
## 
##       Estimate Std. Error z value  Pr(>|z|)    
## ar1  0.0076959  0.3606071  0.0213 0.9829733    
## ma1 -1.2075707  0.3493580 -3.4565 0.0005472 ***
## ma2  0.2075831  0.3485135  0.5956 0.5514262    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

0.4.5.6 ARIMA(2,1,2)

model6.da=Arima(train.diff, order=c(2,1,2),method="ML")
summary(model6.da) #AIC=4749.52

## Series: train.diff 
## ARIMA(2,1,2) 
## 
## Coefficients:
##           ar1      ar2      ma1      ma2
##       -0.9715  -0.2494  -0.2135  -0.7865
## s.e.   0.1406   0.0898   0.1188   0.1179
## 
## sigma^2 = 8707:  log likelihood = -703.32
## AIC=1416.63   AICc=1417.17   BIC=1430.48
## 
## Training set error measures:
##                    ME     RMSE      MAE MPE MAPE      MASE        ACF1
## Training set 1.591334 91.32977 70.49328 Inf  Inf 0.6093032 0.009137627

lmtest::coeftest(model6.da) #terdapat parameter tidak signifikan

## 
## z test of coefficients:
## 
##      Estimate Std. Error z value  Pr(>|z|)    
## ar1 -0.971452   0.140601 -6.9093 4.870e-12 ***
## ar2 -0.249405   0.089781 -2.7779  0.005471 ** 
## ma1 -0.213496   0.118841 -1.7965  0.072418 .  
## ma2 -0.786503   0.117872 -6.6725 2.514e-11 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Berdasarkan pendugaan parameter di atas, nilai AIC terkecil dimiliki oleh model ARIMA(0,1,2) dan parameter model ARIMA(0,1,2) juga seluruhnya signifikan sehingga model yang dipilih adalah model ARIMA(0,1,2).

0.4.6 Analisis Sisaan

Model terbaik hasil identifikasi kemudian dicek asumsi sisaannya. Sisaan model ARIMA harus memenuhi asumsi normalitas, kebebasan sisaan, dan kehomogenan ragam. Diagnostik model dilakukan secara eksplorasi dan uji formal.

0.4.6.1 Eksplorasi Sisaan

#Eksplorasi
sisaan.da <- model4.da$residuals
par(mfrow=c(2,2))
qqnorm(sisaan.da)
qqline(sisaan.da, col = "blue", lwd = 2)
plot(c(1:length(sisaan.da)),sisaan.da)
acf(sisaan.da)
pacf(sisaan.da)

par(mfrow = c(1,1))

Berdasarkan plot kuantil-kuantil normal, secara eksplorasi ditunjukkan sisaan tidak menyebar normal ditandai dengan titik titik yang cenderung tidak mengikuti garis \(45^{\circ}\). Kemudian dapat dilihat juga lebar pita sisaan yang cenderung tidak sama menandakan bahwa sisaan memiliki ragam yang heterogen. Plot ACF dan PACF sisaan ARIMA(0,0,2) juga tidak signifikan pada 20 lag awal yang menandakan saling bebas. Kondisi ini akan diuji lebih lanjut dengan uji formal.

0.4.6.2 Uji Formal

#1) Sisaan Menyebar Normal
ks.test(sisaan.da,"pnorm")  #tak tolak H0 > sisaan menyebar normal

## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  sisaan.da
## D = 0.5042, p-value < 2.2e-16
## alternative hypothesis: two-sided

\(H_0\) : Sisaan menyebar normal

\(H_1\) : Sisaan tidak menyebar normal

Berdasarkan uji KS tersebut, didapat p-value sebesar 0.00 yang kurang dari taraf nyata 5% sehingga tolak \(H_0\) dan menandakan bahwa sisaan tidak menyebar normal. Hal ini sesuai dengan hasil eksplorasi menggunakan plot kuantil-kuantil normal.

#2) Sisaan saling bebas/tidak ada autokorelasi
Box.test(sisaan.da, type = "Ljung")  #tak tolak H0 > sisaan saling bebas

## 
##  Box-Ljung test
## 
## data:  sisaan.da
## X-squared = 8.254e-05, df = 1, p-value = 0.9928

Selanjutnya akan dilakukan uji formal untuk kebebasan sisaan menggunakan uji Ljung-Box. Hipotesis yang digunakan adalah sebagai berikut.

\(H_0\) : Sisaan saling bebas

\(H_1\) : Sisaan tidak tidak saling bebas

Berdasarkan uji Ljung-Box tersebut, didapat p-value sebesar 0.8471 yang lebih besar dari taraf nyata 5% sehingga tak tolak \(H_0\) dan menandakan bahwa sisaan saling bebas. Hal ini berbeda dengan eksplorasi.

#3) Sisaan homogen
Box.test((sisaan.da)^2, type = "Ljung")  #tak tolak H0 > sisaan homogen

## 
##  Box-Ljung test
## 
## data:  (sisaan.da)^2
## X-squared = 0.0018275, df = 1, p-value = 0.9659

Hipotesis yang digunakan untuk uji kehomogenan ragam adalah sebagai berikut.

\(H_0\) : Ragam sisaan homogen

\(H_1\) : Ragam sisaan tidak homogen

Berdasarkan uji Ljung-Box terhadap sisaan kuadrat tersebut, didapat p-value sebesar 0.000 yang kurang dari taraf nyata 5% sehingga tak tolak \(H_0\) dan menandakan bahwa ragam sisaan tidak homogen.

#4) Nilai tengah sisaan sama dengan nol
t.test(sisaan.da, mu = 0, conf.level = 0.95)  #tak tolak h0 > nilai tengah sisaan sama dengan 0

## 
##  One Sample t-test
## 
## data:  sisaan.da
## t = 0.18908, df = 118, p-value = 0.8504
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -15.24941  18.46884
## sample estimates:
## mean of x 
##  1.609717

Terakhir, dengan uji-t, akan dicek apakah nilai tengah sisaan sama dengan nol. Hipotesis yang diujikan sebagai berikut.

\(H_0\) : nilai tengah sisaan sama dengan 0

\(H_1\) : nilai tengah sisaan tidak sama dengan 0

Berdasarkan uji-ttersebut, didapat p-value sebesar 0.4866 yang lebih besar dari taraf nyata 5% sehingga tak tolak \(H_0\) dan menandakan bahwa nilai tengah sisaan sama dengan nol. Hal ini berbeda dengan eksplorasi.

0.4.7 Peramalan

Peramalan dilakukan menggunakan fungsi forecast() . Contoh peramalan berikut ini dilakukan untuk 30 hari ke depan.

#---FORECAST---#
ramalan.da <- forecast::forecast(model4.da, h = 30)
ramalan.da

##     Point Forecast     Lo 80    Hi 80     Lo 95    Hi 95
## 121       9.104968 -111.4603 129.6702 -175.2837 193.4936
## 122     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 123     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 124     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 125     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 126     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 127     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 128     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 129     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 130     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 131     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 132     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 133     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 134     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 135     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 136     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 137     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 138     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 139     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 140     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 141     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 142     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 143     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 144     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 145     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 146     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 147     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 148     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 149     -13.534205 -136.2910 109.2226 -201.2746 174.2061
## 150     -13.534205 -136.2910 109.2226 -201.2746 174.2061

data.ramalan.da <- ramalan.da$mean
plot(ramalan.da)

Berdasarkan hasil plot ramalan di atas, dapat dilihat bahwa ramalan ARIMA(0,012) cenderung stabil hingga akhir periode. Selanjutnya, dapat dicari nilai akurasi antara hasil ramalan dengan data uji sebagai berikut.

pt_1 <- train.ts[430] #nilai akhir data latih
hasil.forc.Diff <- data.ramalan.da
hasil <- diffinv(hasil.forc.Diff, differences = 1) + pt_1
#has.1 sama hasilnta dengan: cumsum(c(pt_1,hasil.forc.Diff))
ts.plot(train.ts,hasil)

perbandingan.da<-matrix(data=c(head(test.ts, n=30), hasil[-1]),
                     nrow = 30, ncol = 2)
colnames(perbandingan.da)<-c("Aktual","Hasil Forecast")
perbandingan.da

##       Aktual Hasil Forecast
##  [1,]   2355             NA
##  [2,]   2356             NA
##  [3,]   2302             NA
##  [4,]   2272             NA
##  [5,]   2324             NA
##  [6,]   2269             NA
##  [7,]   2299             NA
##  [8,]   2284             NA
##  [9,]   2279             NA
## [10,]   2326             NA
## [11,]   2351             NA
## [12,]   2335             NA
## [13,]   2360             NA
## [14,]   2269             NA
## [15,]   2275             NA
## [16,]   2226             NA
## [17,]   2208             NA
## [18,]   2191             NA
## [19,]   2149             NA
## [20,]   2126             NA
## [21,]   2152             NA
## [22,]   2187             NA
## [23,]   2141             NA
## [24,]   2089             NA
## [25,]   2107             NA
## [26,]   2031             NA
## [27,]   1984             NA
## [28,]   1983             NA
## [29,]   1981             NA
## [30,]   2016             NA

accuracy(ts(hasil[-1]), head(test.ts, n=30))

##           ME RMSE MAE MPE MAPE ACF1 Theil's U
## Test set NaN  NaN NaN NaN  NaN   NA         0